1. Field of the Invention
This invention relates to the field of molecular modeling and more particularly to the field of calculating solvation energy.
2. Description to the Related Art
In the computational analysis of molecules, inter-atomic forces are modeled to determine the configurations of individual molecules, multi-molecular entities, and the tendencies of separate molecules to form bound complexes. Because these configuration and affinities are highly dependent on the interaction between the modeled molecules and the surrounding solvent molecules (typically water), a variety of models have been devised to account for the molecule-solvent interactions.
One such method involves explicitly including the individual solvent molecules in the area around the molecules of interest. However, this is very expensive computationally, because thousands of solvent-related atomic positions and inter-atomic forces must often be computed, in addition to those involving the molecules of interest themselves. Accordingly, many attempts have been made to simplify this analysis, and to compute a free energy of solvation for a given configuration of molecule or multi-molecular entity, based on principles of electrostatics that treat the solvent as a dielectric continuum in which the molecules of interest are embedded.
The following formula represents one commonly used model of calculating the solvation energy of a molecule:Gsol=Gpol+Gnp  (Eq.1)
In the above formula, Gsol represents the solvation energy of the molecule, Gpol represents the polarization component of the solvation energy (polar solvation energy), and Gnp represents the non-polar component of the solvation energy (non-polar solvation energy). The non-polar component Gnp can be estimated using the following formula:Gnp=σSA  (Eq.2)
In the above formula, SA is the total solvent accessible area of the molecule, and σ is an empirical parameter. The above equation can be modified to a summation over individual atoms of the molecule, as shown in the following equation:Gnp=ΣiσiSAi  (Eq.3)
The polar solvation energy Gpol can be calculated by methods such as the Poisson-Boltzman (PB) method, the generalized Born (GB) method, the analytical continuum electrostatic potential (ACE) method, and so forth.
In the generalized Born method, the polar solvation energy is determined using the following equation:Gpol=−166(1/∈m−1/∈w)ΣiΣjqiqj/(rij2+aiajexp(−rij2/4aiaj))0.5.  (Eq.4)
In the above formula, ∈m is a dielectric constant of the molecule, ∈w is a dielectric parameter of the solvent, qi is the partial charge of the atom i, qj is the partial charge of the atom j, ai and aj are known as the respective effective Born radii of the atoms i and j, and rij is the distance between atoms i and j. In the above formula, the summations are taken over all atoms of the molecule. If the effective Born radius of every atom is established, then the molecule's polar solvation energy can be calculated.
The effective Born radius of a given atom in a molecule may be generally characterized as the van der Waals radius of the atom increased by an amount characterizing the degree to which the atom is screened from interacting with the solvent by the other atoms of the molecule. More screening equates to a larger effective Born radius. A larger effective Born radius equates to a smaller contribution to the molecule's polar solvation energy computed with Eq.4.
The effective Born radius of an atom i is typically calculated by setting the charge on the atom qi to 1, and approximating the polar component of the solvation energy of the atom Gpol,i of the molecule using a number of numerical methods, or an analytical solution when possible. Once the polar solvation energy for the an atom is estimated, the effective Born radius is calculated as follows:ai=−166(1/∈m−1/∈w)/Gpol,i.  (Eq.5)
Many practical methods of calculating a; and Gpol,i use a Coulomb Field approximation:ai−1=1/4π∫exr−4dV=Ri−1−1/(4π)∫in,r>Rir−4dV  (Eq.6)Gpol,i=−(1/∈m−1/∈w)166/4π∫exr−4dV=−(1/∈m−1/∈w)166(Ri−1−1/4π∫in,r>Rir−4dV)  (Eq.7)
The first term −(1/∈m−1/∈w)166Ri−1 may be termed the “self energy” of the solvated atom if regarded as a single ion having one unit positive charge. The term ∫in,r>Rir−4dV in the above equations Eq.6 and Eq.7 is an integration over the region inside the molecule but outside the van der Waals radius of atom i, which evaluates the screening effect of the rest of the molecule. More details of the above equation are provided in Bashford & Case, Ann. Rev. Phys. Chem., 51, 129–152 (2000), which is incorporated herein by reference in its entirety.
According to a method introduced by Still, the polar component Gpol,i for atom i is calculated using the following discrete sum over the atoms of the molecule rather than a volume integral:Gpol,i=−166(1/∈m−1/∈w)[1/(P0+Ri)−ΣPVj/rij4].  (Eq.8)
In the above equation, Ri is the van der Waals radius of atom i, Vj is the volume of atom j defined by its van der Waals radius, and P0 and P are empirically determined parameters. As shown in the above equations, the polar solvation energy contribution of an atom i is related to its position relative to other atoms j of the molecule.
The determined polar component Gpol,i is used in Eq.6 to determine the effective Born radius for the atom i.
The determined effective Born radii for the atoms are then used in Eq.4 to determine the polar solvation energy of the molecule. With the polar solvation energy determined, the solvation energy can be determined using Eq.1. More details of the generalized Born method and the Still method are provided in U.S. Pat. No. 5,420,805 and Still et al., J. Am. Chem. Soc., 112, 6127–6129 (1990), both of which are hereby incorporated by reference in their entireties.
A variety of modifications and improvements to the original Still formula have been developed which also involve sums over the atoms of the molecule. For example, the following linear pair-wise approximation equation can also be used to derive a polar component of the solvation energy of an atom in a molecule:Gpol,i=−(1/∈m−1/∈w)[166/(λRvdw,i)+166*P1/(R2vdw,i)+Σbond,jP2Vj/rij4+Σangle,jP3Vj/rij4+Σnonbond,jP4Vj/rij4C(rij, Rvdw,i, Rvdw,j, P5)].  (Eq.9)
In the above equation, λ, P1, P2, P3, P4, and P5 represent empirical parameters. Rvdw,i and Rvdw,j represent the van der Waals radii of atoms i and j. The summation terms are used in Eq.9 to approximate the polar component of an atom i. In Eq.9, the other atoms j are separated into three types: atoms that form direct bonds with atom i (represented by the Σbond,j term), atoms that form indirect bonds with atom i (represented by the Σangle,j term), and atoms that do not form bonds with atom i (represented by the Σnonbond,j term). C(rij, Rvdw,i, Rvdw,j, P5) represents an analytical function that is at least a function of rij, Rvdw,i Rvdw,j, and P5. The other atoms j include all atoms of the molecule that contribute to the polar component of the atom i. In one embodiment for ease of formulation, the other atoms j include all atoms of the molecule. More details of the above equation are provided in Dominy & Brooks, J. Phys. Chem. B, 103, 3765–3773 (1999), incorporated herein by reference in its entirety.
Another pair-wise approximation equation similar to Eq.9 can also be used:Gpol,i=−(1/∈m−1/∈w)[166/(Rvdw,I+φ+P1)+Σbond,jP2Vj/rij4+Σangle,jP3Vj/rij4+Σnonbond,jP4Vj/rij4C(rij, Rvdw,i, Rvdw,j, P5)].  (Eq.10)
In the above equation, φ, P1, P2, P3, P4, and P5 represent empirical parameters. C(rij, Rvdw,i, Rvdw,j, P5) represents an analytical function that is at least a function of rij, Rvdw,i, Rvdw,j, and P5. More details of the above equation are provided in Qui et al., J. Phys. Chem. B., 101, 3005–3014 (1997), also incorporated herein by reference in its entirety. As described above with respect to the integral formulations, the first term or terms of these equations which do not involve sums over the atoms of the system can be considered the “self energy” terms, and the sums are the “screening effect” terms.
Compared to explicitly including solvent molecules in the simulation, the generalized Born method requires fewer calculations, while producing results similar to those obtained by a more computationally intensive method, such as by using the Poisson-Boltzmann equation. Therefore, the generalized Born method facilitates dynamic simulation with long trajectories for large molecule systems, such as globular proteins and DNA.
However, the existing generalized Born method assumes that the solvent is homogeneous in all directions surrounding the molecular system being modeled. It has thus not been applied to study solvation effects in non-homogeneous media, or in environments having regions outside the solute molecule(s) which have different dielectric constants. As mentioned above, one example of such a system with biological and clinical significance are proteins which are partially embedded in cell membrane and partially in aqueous solution. Trans-membrane proteins are important drug targets, but since trans-membrane proteins are typically partially embedded in a hydrophobic lipid bilayer of cell membrane and partially embedded in a polar aqueous solution, the generalized Born methods described above have not been suitable for modeling trans-membrane proteins. Furthermore, since methods of crystallizing trans-membrane proteins are under-developed, and the results from solid state NMR measurements are still sparse, little is known about the structure of these proteins in native form and the mechanism of their actions. The limited number of 3-D structures limits the use of the homology modeling method for predicting trans-membrane protein structures.
Although other computational methods such as the Poisson-Boltzmann equation can be used to calculate the solvation energy of a molecule in a more complex solvent system, the computational complexity is significant.